Hepatocyte population analysis

General Setup

Setup chunk

Setup reticulate

knitr::opts_chunk$set(fig.width = 8)
knitr::opts_knit$set(root.dir = normalizePath(".."))
knitr::opts_knit$get("root.dir")

[1] "/nas/groups/treutlein/USERS/tomasgomes/projects/liver_regen"

Load libraries

library(reticulate)
knitr::knit_engines$set(python = reticulate::eng_python)
py_available(initialize = FALSE)

[1] FALSE

use_python(Sys.which("python"))
py_config()

python:         /home/tpires/bin/miniconda3/bin/python
libpython:      /home/tpires/bin/miniconda3/lib/libpython3.8.so
pythonhome:     /home/tpires/bin/miniconda3:/home/tpires/bin/miniconda3
version:        3.8.3 (default, May 19 2020, 18:47:26)  [GCC 7.3.0]
numpy:          /home/tpires/bin/miniconda3/lib/python3.8/site-packages/numpy
numpy_version:  1.18.5

NOTE: Python version was forced by RETICULATE_PYTHON

Load data (from all cells)

library(Seurat)

Attaching SeuratObject

library(ggplot2)
library(destiny)
library(plyr)
library(dplyr)

Attaching package: ‘dplyr’

The following objects are masked from ‘package:plyr’:

    arrange, count, desc, failwith, id, mutate, rename, summarise, summarize

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

library(RColorBrewer)

Immune Hepatocytes

Subset Hepatocytes

allcells_css = readRDS(file = "data/processed/allcells_css.RDS")

Choose dimensionality reduction

# the dividing cells are T/NK, will update label at the end
hep_pops = unique(allcells_css$allcells_clusters)[grepl("Hepatocytes ",
                                                        unique(allcells_css$allcells_clusters))]
all_hep_cells = allcells_css[,allcells_css@meta.data$allcells_clusters %in% hep_pops]
all_hep_cells = suppressWarnings(SCTransform(all_hep_cells, do.correct.umi = T, verbose = F, 
                                             vars.to.regress=c("unique_name","nCount_RNA"),
                                             variable.features.rv.th = 1, seed.use = 1,
                                             return.only.var.genes = F, 
                                             variable.features.n = NULL))

all_hep_cells = RunPCA(all_hep_cells, verbose = F)

Cluster and get markers

pcadonorplt = list()
for(i in 1:15){
  n = 2*i
  pcadonorplt[[i]] = DimPlot(all_hep_cells, dims = (n-1):n, 
                             reduction = "pca", group.by = "Donor")+
    theme(aspect.ratio = 1, legend.position = "none")
}
cowplot::plot_grid(plotlist = pcadonorplt, ncol = 5, align = "hv")

ldims = list(1:30, c(1:2,4:30), c(1:2,4:20), c(1:2,4:16), c(2,4:16, 18, 20, 24:26))
for(i in 1:length(ldims)){
  all_hep_cells = RunUMAP(all_hep_cells, dims = ldims[[i]], verbose = F)
  plt1 = DimPlot(all_hep_cells, reduction = "umap", group.by = "Condition")
  plt2 = DimPlot(all_hep_cells, reduction = "umap", group.by = "Donor")
  print(plt1+plt2)
}

all_hep_cells = RunUMAP(all_hep_cells, dims = c(2,4:16, 18, 20, 24:26), verbose = F)
DimPlot(all_hep_cells, reduction = "umap", group.by = "Condition")

DimPlot(all_hep_cells, reduction = "umap", group.by = "Donor")

DimPlot(all_hep_cells, reduction = "umap", group.by = "Phase")

FeaturePlot(all_hep_cells, reduction = "umap", ncol = 4,
            features = c("ALB", "SAA1", "SAA2", "HAMP", "CYP1A2", "CYP2E1", "NNMT", "MT1G"))

Save data

DimPlot(all_hep_cells, reduction = "umap", group.by = "pca20_res.0.8", label = T)
DimPlot(all_hep_cells, reduction = "umap", group.by = "pca20_res.0.3", label = T)

DimPlot(all_hep_cells, reduction = "umap", group.by = "Donor", label = F)

DimPlot(all_hep_cells, reduction = "umap", group.by = "Condition", label = F)

annot_fine = list("10" = "Hepatic cancer/regeneration response", # https://journals.sagepub.com/doi/pdf/10.1177/153537020623100203, https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2292801/, https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4656314/, https://onlinelibrary.wiley.com/doi/full/10.1002/cbf.3288
                  "12" = "Dividing cells",
                  "11" = "Hepatocyte-Immune interaction",
                  "7" = "Periportal inflammation-inducing", # https://www.sciencedirect.com/science/article/pii/S0753332220308738, 
                  "0" = "Hepatocytes (HAMP+)",
                  "1" = "Pericentral Hepatocytes",
                  "8" = "Pericentral regenerating", # https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5349637/
                  "3" = "Fibrotic and regenerating response", # https://pubmed.ncbi.nlm.nih.gov/23770341/, https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5937820/
                  "5" = "Unknown 1", # https://onlinelibrary.wiley.com/doi/full/10.1002/cbf.3288
                  "2" = "Ribosomal",
                  "4" = "Unknown 2"
                  )

round(t(t(table(all_hep_cells$pca20_res.0.8, all_hep_cells$Condition))/colSums(table(all_hep_cells$pca20_res.0.8, all_hep_cells$Condition)))*100, 2)

round(t(t(table(all_hep_cells$pca20_res.0.3, all_hep_cells$Condition))/colSums(table(all_hep_cells$pca20_res.0.3, all_hep_cells$Condition)))*100, 2)

round(t(t(table(all_hep_cells$pca20_res.0.3, all_hep_cells$Donor))/colSums(table(all_hep_cells$pca20_res.0.3, all_hep_cells$Donor)))*100, 2)

pheatmap::pheatmap(round(t(t(table(all_hep_cells$pca20_res.0.8, all_hep_cells$Condition))/colSums(table(all_hep_cells$pca20_res.0.8, all_hep_cells$Condition)))*100, 2))
pheatmap::pheatmap(round(t(t(table(all_hep_cells$pca20_res.0.3, all_hep_cells$Condition))/colSums(table(all_hep_cells$pca20_res.0.3, all_hep_cells$Condition)))*100, 2))